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ABSTRACT 

We show that the globular cluster mass function (GCMF) in the Milky Way depends on cluster half-mass 
density, p/,, in the sense that the turnover mass Mjo increases with p/, while the width of the GCMF decreases. 
We argue that this is the expected signature of the slow erosion of a mass function that initially rose towards low 
masses, predominantly through cluster evaporation driven by internal two-body relaxation. We find excellent 
agreement between the observed GCMF — including its dependence on internal density r/zo/,, central concen- 
tration c, and Galactocentric distance r„^ — and a simple model in which the relaxation-driven mass-loss rates 

1 /2 

of clusters are approximated by -dMjdt = /iev oc . In particular, we recover the well-known insensitivity 
of Mjo to rgc. This feature does not derive from a literal "universality" of the GCMF turnover mass, but rather 
from a significant variation of Mxo with pi, — the expected outcome of relaxation-driven cluster disruption — 
plus significant scatter in p/, as a function of r^^. Our conclusions are the same if the evaporation rates are 
assumed to depend instead on the mean volume or surface densities of clusters inside their tidal radii, as 

Mev oc p^^ or /iev Y}^^^ — alternative prescriptions that are physically motivated but involve cluster properties 
(p, and E,) that are not as well defined or as readily observable as p/,. In all cases, the normalization of p^y 
required to fit the GCMF implies cluster lifetimes that are within the range of standard values (although falling 
towards the low end of this range). Our analysis does not depend on any assumptions or information about 
velocity anisotropy in the globular cluster system. 

Subject headings: galaxies: star clusters — globular clusters: general 



1. INTRODUCTION 

The mass functions of star cluster systems provide an im- 
portant point of reference for attempts to understand the con- 
nection between old globular clusters (GCs) and the young 
massive clusters that form in local starbursts and galaxy merg- 
ers. When expressed as the number per unit logarithmic mass, 
dN/d log M, the GC mass function (GCMF) is character- 
ized by a peak, or turnover, at a mass Mxo ~ 1-2 x 10^ Mq 
that is empirically very similar in most galaxies. By con- 
trast, the mass functions of young clusters show no such fea- 
ture but instead rise monotonically towards low masses over 
the full observed range (10'' Mq >M> IO'^Mq in the best- 
studied cases), in a way that is well de scribed by a power l aw, 
dN/d log M oc M'-^ with /3 ~ 2 (e.g. JZhang & Falllll999h . 

At the same time, for high M > Mjq, old GCMFs 
closely resemble the mass functions of young clusters, and 
of molecular clouds in the Milky Way and other ga laxies 
(iHarris & Pudritd Il994t lElmegreen & Efremovl 1 19971) : and 
it is well known that a number of dynamical processes 
cause star clusters to lose mass and can lead to their com- 
plete destruction as they orbit for a Hubble time in the 
potential wells of their parent galaxies (e.g.. Fall & Rees 
1977; Caputo & Castellani 1984; Aguilar, Hut, & Ostriker 
1988; ,Chernoff& Weinberg., 1990: Gnedin & Ostriker .1997: 



Murali & WeinbergI 1 19971) . It is therefore natural to ask 



whether the peaks in GCMFs can be explained by the deple- 
tion over many Gyr of globulars from initial mass distribu- 
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tions that were similar to those of young clusters below Mxo 
as well as above. 

Our chief purpose in this paper is to establish and inter- 
pret an aspect of the Galactic GCMF that appears fundamen- 
tal but has gone largely unnoticed to date: dN/d log M has 
a strong and systematic dependence on GC half-mass den- 
sity. Pi, = 3M /%nr\ (r/, being the cluster half-mass radius), in 
the sense that the turnover mass Mxo increases and the width 
of the distribution decreases with increasing pi,. As observed 
facts, these must be explained by any theory of the GCMF. We 
argue here that they are an expected signature of slow dynam- 
ical evolution from a mass function that initially increased 
towards M < Mxo, if the long-term mass loss from surviv- 
ing GCs has been dominated by stellar escape due to internal, 
two-body relaxation (which we refer to from now on as either 
re laxation-driven evapo ration or simply evaporation). 

iFall & Zhang! (1200 H hereafter IfZOiI) explain in detail why 
cluster evaporation dominates the long-term evolution of the 
low-mass shape of observable GCMFs. Briefly, stellar evolu- 
tion removes (through supernovae and winds) the same frac- 
tion of mass from all clusters of a given age, and so cannot 
change the shape of d N/d log M (unless speci al initial con- 
ditions are invoked; cf. IVesperini & Zepll2003 '). Meanwhile, 
for GCs like those that have survived for a Hubble time in 
the Milky Way, the mass loss from gravitational shocks dur- 
ing disk crossings and bulge passag es is generally less than 
that d ue to evaporation for M < Mxo (iFZOlbl^eto & Gnediiil 
l2006l).5 

As we discuss further in ^below, the evaporation of tidally 
limited clusters proceeds at a rate, p^^ = -dM/dt, that is ap- 
proximately constant in time and primarily determined by 

^ It is possible that there existed a past population of GCs with low den- 
sities or concentrations, or perhaps on extreme orbits, that were destroyed in 
less than a Hubble time by shocks or stellar evolution. Our discussion does 
not cover such clusters. 
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cluster density. iFZOll show that a constant mass-loss rate 
leads to a power-law scaling dN /d log M oc M'^~^ with /3 — ^ 
(corresponding to a flat distribution of clusters per unit lin- 
ear mass) at sufficiently low M < ji^^t in the evolved mass 
function of coeval GCs that began with any nontrivial ini- 
tial dN /d log Mq.^ To accommodate this when dN /d log Mq 
originally increased towards low masses as a power law, a 
time-dependent peak must develop in the GCMF at a mass 
of order Mjo(0 ~ Mev? (FZOl). But then, since /iev depends 
fundamentally on cluster density, so too must Mto- 

A /3 ~ power-law scaling below the turnover mass has 
been c onfirmed directly in the G CMFs of the Milk y Way 
(iFZOlh and the gian t elliptical M87 (IWaters et al.l2006l) . while 
[Jordan et al] (l2007h show it to be consistent with dN /d log M 
data for 89 Virgo Cluster galaxies, and it is apparent in deep 
observations o f some other GCM Fs (e.g., in the Sombrero 
galaxy, M104: ISpitler et al.ll2006h. As regards the peak itself, 
old GCs are observed (e.g. . [Jordan et al.ll2005l) to have rather 
similar densities on average — and, therefore, similar typical 
/iev — in galaxies with widely different total luminosities and 
Hubble types. (Inasmuch as cluster densities are set by tides, 
this is probably related to the mild variation of mean galaxy 
density with total luminosity; see FZOl, and also Jordan et al.l 
l2007i ) Thus, an evaporation-dominated evolutionary origin 
for a turnover in the GCMF appears to be consistent with 
the well-known fact that the mass scale Mtq generally dif- 
fers v ery little among galaxies (e.g.. lHarrisll2001l: [Jordan et al.l 
l2006h . 

If this picture is basically correct, it implies that, even 
though Mto may appear nearly universal when considering 
the global mass functions of entire GC systems, in fact the 
GCMFs of subsamples of clusters with similar ages but dif- 
ferent densities should have different turnovers. In ^ we 
show — working for definiteness and relatively easy observ- 
ability with the half-mass density, pi, — that this is the case for 
globulars in the Milky Way. We fit the observed dN /d log M 
for GCs in bins of different /?/, with models assuming that 
(1) the initial distribution increased as a /3 = 2 power law at 
low masses and (2) the mass-loss rates of individual clusters 
can be estimated from their half-mass densities by the rule 

Mev p]!'^ . In ^we discuss the validity of this prescription 
for /iev, which is certainly approximate but captures the main 
physical dependence of relaxation-driven mass loss. In par- 

1 /2 

ticular, we show that the alternative mass-loss laws /iev p, ' 

and /iev Y^J'^ — where p, and S, are the mean volume and 
surface densities inside cluster tidal radii — lead to models for 
the GCMF that are essentially indistinguishable from those 

1/2 

based on /iev « Pk ■ The normalization of /igv required to fit 
the observed GCMF implies cluster lifetimes that are within 
a factor of « 2 (perhaps slightly on the low side, if the ini- 
tial power-law exponent at low masses was /3 = 2) of typical 
values in theories and simulations of two-body relaxation in 
tidally limited GCs. 

We also show in ^that when the observed densities of in- 
dividual clusters are used in our models to predict GCMFs 
in different bins of Galactocentric radius (r^c), they fit the 
much weaker variation of dN/d log M and Mto as functions 
of Tgc, which is well-known in the Milky Way and other large 

* Throughout this paper, we use "initial" to mean at a relatively early time 
in the development of long-lived clusters, after they have dispersed any rem- 
nants of their natal gas clouds, survived the bulk of stellar-evolution mass 
loss, and come into virial equilibrium in the tidal field of a galaxy. 



galaxies (s ee iHarrisI 120011: [Harris. Har ris. & McLaughhnI 
1998; Barm bv. Huchra. & Brodiel2001l: IVesperini et al.ll2003r 
iJordan et al.l l2007h . Similarly, applying our models to the 
GCs in two bins of central concentration, with only the mea- 
sured ph of the clusters in each subsample as input, suffices 
to account for previously noted differences between the mass 
functions of low- and h igh-concentration Galactic globulars 
jSmith & Burkertl2002l) . The most fundamental feature of the 
GCMF therefore appears to be its dependence on cluster den- 
sity, which can be understood at least qualitatively (and even 
quantitatively, to within a factor of 2) in terms of evaporation- 
dominated cluster disruption. 

There is a widespread perception that if the GCMF evolved 
slowly from a rising power law at low masses, then a weak 
or null variation of Mto with rgc can be achieved only in GC 
systems with strongly radially anisotropic v elocity distribu- 
tions, which are not observed (see especially IVesperini et al.l 
l2003h . This apparent inconsistency has been cited to bol- 
ster some recent attempts to identify a mechanism by which 
a "universal" peak at Mto ^ 10^ M0 might have been im- 
printed on the GCMF at the time of cluster formation, or 
very shortly afterwards, and little af fected by the subsequent 
destruction of lower-mass GCs (e.g.. IVesperini & Zepfll200l 
iParmentier & Gilmorell2007h . However, given the real suc- 
cesses of an evaporation-dominated evolutionary scenario for 
the origin of Mto, as summarized above and added to below, 
it would be premature to reject the idea in favor of requiring a 
near-formation origin, solely on the basis of difficulties with 
GC kinematics. (And, in any event, formation-oriented mod- 
els must now be reconsidered in light of the non-universality 
of Mto as a function of cluster density.) 

We are not concerned in this paper with velocity anisotropy 
in GC systems, because we only predict an evaporation- 
evolved dN/d log M as a function of cluster density (and age) 
and take the observed distribution of p/, versus rgc in the Milky 
Way as a given, to show consistency with the observed be- 
havior of Mto as a fu nction of rgc. Most other models (FZOU 
IVesperini et al.ll2003t and references therein) predict dynam- 
ically evolved GCMFs directly in terms of r^c, and in doing 
so are forced also to derive theoretical dependences of cluster 
density on rgc. It is only at this stage that GC orbital distri- 
butions enter the problem, and then only in conjunction with 
several other assumptions and simplifications. As we discuss 
further in |3]b6low, the radially biased GC velocity distribu- 
tions that appear in such models could well be consequences 
of one or more of these other assumptions, rather than of the 
main hypothesis about evaporation-dominated GCMF evolu- 
tion. 

2. THE GALACTIC GCMF AS A FUNCTION OF 
CLUSTER DENSITY 

In this section we define and model the dependence of 
the Galactic GCMF on cluster density. First, we describe 
the dependence that is expected to arise from evaporation- 
dominated evolution. 

Two-body relaxation in a tidally limited GC leads to a 
roughly steady rate of mass loss, /tev = -dM/dt ~ constant 
in time. Thus, the total cluster mass decreases approxi- 
mately linearly, as M(f) ~ Mq- p^yt. Thi s behavior is exact 
in some classic models of GC evolution (lHenon|[l961l) and 
is found t o be a good approximation in most other calcula- 
tions (e.g.,|L ee & Ostriker 1987; Chernoff & Weinberg 199C; 
Vesperini & Heggie 1997; Gnedin. Lee. & Osti'iker 19991; 



Baumgardtll200iriGierszll2001l: iBaumgardt & Makinoll2003 
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iTrenti. Heggie. & Hutl2007h . The result comes from a variety 
of computational methods (semi-analytical, Fokker-Planck, 
Monte Carlo, and A^-body simulation) applied to clusters with 
different initial conditions (densities and concentrations) on 
different kinds of orbits (circular and eccentric; with and with- 
out external gravitational shocks) and with different internal 
processes and ingredients (with or without stellar mass spec- 
tra, binaries, and central black holes). To be sure, deviations 
from perfect linearity in M(t) do occur, but these are generally 
small — especially away from the endpoints of the evolution, 
i.e., forO.9 > M(f)/Mo > 0.1 — and neglecting them to assume 
an approximately constant dM/dt is entirely appropriate for 
our purposes. 

When gravitational shocks are subdominant to relaxation- 
driven evaporation, as they generally appear to be for 
extant GCs, they work to boost the mass-loss rate 
sUghtly without altering the basic Unearity of M(t) ( e.g., 
IVesperini & Hegpelll997l:lGnedin. Lee. & OstrikeHI 19991: see 
also Figure 1 of FZOlh . A time-dependent mass scale A = 
/ievf is then associated naturally with any system of coeval 
clusters having a common mass-loss rate: all those with initial 
Mo < A are disrupted by time t, and replaced with the rem- 
nants of objects that began with Mq > A. As we mentioned 
in ^ if the initial GCMF increased towards low masses as a 
power law, then A is closely related to a peak in the evolved 
distribution, which eventually decr eases t owards low M < A 
as dN/d log M oc M'"^ wi th /3 = (iFZOll). 

In s tandard theory (e.g.. lSpitzej|1987l : [Binney & Tremaind 
119871 Section 8.3), the lifetime of a cluster against evapora- 
tion is a multiple of its two-body relaxation time, fi-ix. For 
a total mass M of stars within a radius r, this scales to first 
order (ignoring a weak mass dependence in the Coulomb log- 
arithm) as ti-i^(r) oc (Mpy^^ ocM/p'/^, where poiM/P. In a 
concentrated cluster with an internal density gradient, tyi^{r) 
of course varies throughout the cluster, and the global re- 
laxation timescale i s an averag e of the local values (see the 
early discussion by lKing|ll958l) . This can still be written as 
fiix oc M/p'/^, with M the total cluster mass and p an appro- 
priate reference density. We then have for the instantaneous 
mass-loss rate, = —dM/dt ccM/t^\^ (x'p^l^. Insofar as this 
is approximately constant in time, a GCMF evolving from an 
initial /3 > 1 power law at low masses should therefore de- 
velop a peak at a mass that depends on cluster density and age 
through the parameter A oc p'''^?. 

It remains to identify the best measure of p in this context. 
A standard choice in the literature, and the one that we even- 
tually make to derive our main results in this paper, is the 
half-mass density pi, = 3M/87rr,^. However, in a steady tidal 
field, the mean density p, inside the tidal radius of a cluster is 
constant by definition, and thus choosing p = p, instead is the 
simplest way to ensure that /iev oc p'^^ and /iev — constant in 
time are mutually consistent. In fact. King ( 1966) found from 
direct calculations of the escape rate at each radius within his 
standard (lowered Maxwellian) models, that the coefficient 

1 12 

in p, is only a weak function of the internal density 

structure (concentration) of the models, and thus only a weak 
function of time for a cluster evolving quasistatically through 

a series of such models. 

1 11 

The rule pev oc p, is routinely used to set the GC mass- 
loss rates in models for the dynamical evolution of the GCMF, 
although such studies normally express pev immediately in 
terms of orbital pericenters, rp, most often by assuming p, oc 
as for GCs in galaxies whose total mass distributions fol- 



low a singular isothermal sphere ( e-g-, IVesperinilll997l 1 19981 
l2(M I2OQII fVesperini et alJlMia feaumgardt 1998: FZOl). 
This bypasses any explicit examination of the GCMF as a 
function of cluster density, which is our main goal in this pa- 
per. But it is done in part because tidal radii are the most 
poorly constrained of all structural parameters for GCs in the 
Milky Way (their theoretical definition is imprecise and their 
empirical estimation is highly model-dependent and sensitive 
to low-surface brightness data), and they are exceedingly dif- 
ficult if not impossible to measure in distant galaxies. We 
deal with this here by focusing on the GCMF as a function of 
cluster density pi, inside the less ambiguous, empirically bet- 
ter determined and more robust half-mass radius, asking how 

1/2 

simple models with pev oc p^ fare against the data. 

1/2 1/2 
Taking pev oc p^ in place of pev oc p, , which we do to 

construct evaporation-evolved model GCMFs in ^2.21 is most 
appropriate if the ratio p, / p/, is the s ame for all clusters and 
constant in time. This is the case in iHenonI ' s 1^6 1) model 
of GC evolution, and in this limit (adopted by IfZOI. in their 
models for the Galactic GCMF) our analysis is rigorously jus- 
tified. However, real clusters are not homologous {pij Ph dif- 
fers among clusters) and they do not evolve self-similarly (p;, 
may vary in time even if p, does not). The key assumption in 
our models is that pev is approximately independent of time 
for any GC, which is well-founded in any case. By using cur- 
rent Ph values to estimate pev, we do not suppose that the half- 
mass densities are also constant, but we in effect use a single 
number for all GCs to represent a range of (p,/ p/,)'''^- Equiv- 
alently, we ignore a dependence on cluster concentration in 

the normalization of pev c)c p^^. As we discuss further in ^ 
it is reasonable to neglect this complication in a first approx- 
imation because (pr/p/,)'''^ varies much less among Galactic 
globulars than p, and p/, do separately. We demonstrate this 
explicitly by repeating our analysis with p/, replaced by p, and 
recover essentially the same results for the GCMF. 

In ^ we also discuss some recent results, which indicate 
that the timescale for relaxation-driven ev aporatio n depends 
on a slightly less-than-line ar power of frix dBaumg ardt 2001]; 
Baum gardt & Makinoll200 3). We point out that this implies 
that Pev may increase as a modest power of the average sur- 
face density of a cluster as well as (or, in an important special 
case, instead of) the usual volume density. However, we show 
in detail that making the appropriate changes throughout the 
rest of the present section to reflect this possibility does not 
change any of our conclusions. 

2.1. Data 

Figure [T] shows the distribution of mass against half-mass 
density and against Gala ctocent r ic radi us for 146 Milky Way 
GCs in the catalogue of iHarrisI (Il996h .^ along with the dis- 
tribution of Ph versus rgc linking the two mass plots. The 
iHarrisI catalogue actually records the absolute V magni- 
tudes of the GCs. We obtain masses from these by apply- 
ing the population-s ynthesis model mass-to-light ratios 
computed by iMcLaughlin & van der Marell (l2005h for indi- 
vidual clusters based on their metallicities and an assumed 
age of 1 3 Gyr. However, we first multiplied all of the 
iMcLaughhn & van der Marell Ty values by a factor of 0.8 so 
as to obtain a median Ti/ ~ 1 .SM© Lq in the end,*^ consistent 

' Feb. 2003 version; see http://physwww.mcmaster.ca/~har ris/mwgc.dat] . 
Throughout this paper, we use xto denote the median of any quantity x. 
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Fig. 1. — Left: Mass versus three-dimensional half-mass density, pj, = SM/Sirr^^, and versus Galactocentric radius, r^c, for 146 Milky Way GCs in the 

catalogue of lHarrisI )199(i) . The dashed line in the first panel is A/ cx: p,'/^, a locus of approximately constant lifetime against evaporation. Right: Half-mass 
density versus rgc for the same clusters. 



with direct dynamical est imates (see lMcLaughlin|| 2000l an d 
iMcLaughlin & van derM arel 2005; also Barmby et al.l2007h . 

By assigning mass-to-light ratios to GCs in this way, we 
allow for expected differences between clusters with differ- 
ent metallicities. Our application of a corrective factor to 
the population-synthesis values, T^°'', is motivated empiri- 
cally by the fact that their distribution among Galactic GCs 
is strongly peaked around a median Ty°'' ~ 1 .9 Mq Lq , while 
the observed (dynamical) Ty^" lie in a fairly narrow range 

aroun d fy^'^ ~ 1.5 MqLq' dMcLaughlin & van der 
l2005h . However, it is worth noting that the size of this differ- 
ence is similar to what is found in some numerical simulations 
of two-body relaxation over a Hubble time in clusters with a 
spectrum of stellar masses (e.g..'Baumgardt & Makinol2003h . 
In such simulations, T'^" falls below T^°'' due to the prefer- 
ential escape of low-mass stars with high individual M^,/L^, 
(population-synthesis models do not incorporate this or any 
other stellar-dynamical effect). Thus, a median T'^" < T^°'' 
may itself be a signature of cluster evaporation. We might 
then also expect that more dynamically evolved clusters — that 
is, those with shorter relaxation times — could have systemat- 
ically lower ratios of T'^" /T'^'^. However, this is a relatively 
small effect, which is not well quantified theoretically and 
is not clearly evident in real data (t he numbers published by 
IMcLaughlin & van der Marell2005l show no significant corre- 
lation between 'Y'^y^YW ^nd frh for Galactic globulars). We 
therefore proceed, as stated, with a single y^^"/T^°'' = 0.8 
assumed for all GCs. 

Harris (1996) gives the projected half-light radius /?/, for 
141 of the clusters with a mass estimated in this way, and for 
these we obtain the three-dim ensional half-m ass radius from 
the general rule r/, = (4/3)/?/, (ISpitzerlll987h . which assumes 
no internal mass segregation. The remaining five objects have 
mass estimates but no size measurements. To each of these 
clusters, we assign an r/, equal to the median value for those 
of the other 141 GCs having masses within a factor two of the 
one with unknown r/,. In all cases, the half-mass density is 
p,, = 3M/8Trrl. 

The leftmost panel in Figure [T] shows immediately that 
the cluster mass distribution has a strong dependence on 
half-mass density: the median M increases with p/, while 



the scatter in logM — that is, the width of the GCMF — 
decreases. The first of these poi nts is related to the fact that 
r/, correlates poorly withM (e.g.. lDiorgovski & Mevlanll994t 
McLaughlin 2000). The second point, that the dispersion of 
dN/d log M decreases \yith increasing p;,, is behind the find- 
ing (Kavelaars & Hanesll 19971; lGnedinll l997h that the GCMF 
is b road er at very large Galactocentric radii. We return to this 
in^ 

A natural concern, when plotting M against pi, as we have 
done here, is that any apparent correlation might only be a 
trivial reflection of the definition pi, (xM /r\. This may seem 
particularly worrisome because, as we just mentioned, it is 
known that size does not correlate especially well with mass 
for GCs in the Milky Way (or, indeed, in other galaxies). 
However, the lack of a tight M-r/, correlation does not imply 
that all GCs have the same r/,, even within the unavoidable 
measurement errors. The root-mean-square (rms) scatter of 
log r/, about its average value is ±0.3 for Galactic GCs, and 
the 68-percentile spread in log r/, is slightly greater than 0.5, 
or more than a factor of 3 in linear terms (from the data in 
Harris 1996; see, e.g.. Figure 8 of McLaughlin 2000). This 
compares to an rms random measurement error (from formal, 
fitting uncertainties) of (5(log r/,) « 0.05, or about 10% rel- 
ative error; and an rms systematic measurement error (i.e., 
differences in the r/, inferred from fitting different structural 
models to a single cluster) of perha ps 6{logri,) < 0.03; see 
IMcLaughlin & van der Marell ( l2005h . Most of the scatter in 
plots of observed half-light radius versus mass is therefore 
real and contains physical information. The left-hand panel 
of Figure [T]displays this information in a form that highlights 
clear, nontrivial overall trends requiring physical explanation. 

The dashed line in the plot of mass against density traces 

1 /2 T 

the proportionality M cx p,/ , or Mrf^ = constant. Insofar as 
the half-mass relaxation time scales as frh oc (Mr^^)'/^, and 

1 /2 

to the extent that /iev Af/^rh cx pj^' approximates the av- 
erage rate of relaxation-driven mass loss, this line is one of 
equal evaporation time. That such a locus nicely bounds 
the lower envelope of the observed cluster distribution is 
itself a strong hint that relaxation-driven cluster disruption 
has significantly modified the GCMF at low masses (re- 
call that Mrl = constant defines one side of the GC "sur- 
vival triangle" when the M-pi, plot is recast as r;, versus M: 
iFaU & Reesl[T977l: lOkazaki & TosalfT99l lOstriker & Cmediiil 
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TABLE 1 

Milky Way GC Properties in Bins of Density and Galactocentric Radius 



Bin N PI,'-" r,c" M„,in M^ax M" Mjo^ 

[M0PC-3] [kpc] [Mq] [Mq] [Mq] [Mq] 

pli bins 

0.034 < Pa < 76.5 Mq pc-3 48 8.48 12.9 5.63 x 10^ 8.84 x 10' 4.12x10* 3.98x10'' 

78.8 < PA < 526 Mq 49 232 5.6 8.37 x 10^ 1.67x10'' 1.22 x lO' 1.58 xlO'^ 

579 < PA < 5.65 x 10" Mq pc"^ 49 973 3.2 1.93 x 10** 1.30 x 10*" 2.82 x lO"* 2.88 x 10^ 

r„. bins 



0.6 < rgc < 3.2 kpc 47 597 1.9 4.47 x 10' 1.02x10*^ 1.15 x 10^ 2.14x10^ 

3.3<rgc<9.4kpc 50 261 5.2 2.02 x 10^ 1.67x10'^ 1.27 x 10^ 1.66x10^ 

9.6 < rgc < 123 kpc 49 18.4 18.3 5.63 x 10^ 1.30 xlO'' 7.42 x 10"* 8.71x10" 

'■^ The notation £ represents the median of quantity x}' Mto is the peak mass of the model GCMFs traced 
by the soUd curves in each panel of Figure (2] which are given by equation ([3} of the text with (3 = 2, 
Mc =10^ Mq, and individual A given by the observed ph of each cluster through equation 



119971: iGnedin & Ostrikeii [19971) . It is also further evidence 
that the weak correlation of observed r/, with M is due to sig- 
nificant and real differences in cluster radii, since if r/, were 
intrinsically the same for all GCs, then we would see M oc ph 
instead. 

The middle panel of Figure [T] shows the well-known result 
that the typical GC mass depends weakly if at all on Galacto- 
centric radius, at least until large > 30-40 kpc, where there 
are too few clusters to discern any trend. The right-hand panel 
of the figure shows why this is true even though the GCMF 
depends significantly on cluster density: although there is a 
coiTelation between half-mass density and Galactocentric po- 
sition, the large scatter about it is such that convolving the 
observed M versus ph with the observed ph versus rgc results 
in an almost null dependence of M on rgc. 

We now divide the GC sample in Figure [T] roughly into 
thirds, in two different ways: first on the basis of half-mass 
density, and second by Galactocentric radius. These p/, and 
rgc bins are defined in Table [T] which also gives a few sum- 
mary statistics for the globulars in each subsample. We count 
the clusters in every subsample in about 10 equal-width bins 
of log M to obtain histogram representations of dN /d log M, 
first as a function of p/, and then as a function of rgc. These 
GCMFs are shown by the points in Figure |2] with errorbars 
indicating standard Poisson uncertainties. The curves in the 
figure trace model GCMFs, which we describe in ^2.21 For 
the moment, it is important to note that the dashed curve is the 
same in every panel, apart from minor differences in normal- 
ization, and is proportional to the GCMF for the whole sample 
of 146 GCs. (In the middle-left panel of Figure|2l which per- 
tains to clusters distributed tightly around the median pk of 
the entire GC system, the dashed curve is coincident with the 
solid curve running through the data.) 

The left-hand panels of Figure |2] show directly that the 
GCMF is peaked for clusters at any density, and that the mass 
of the peak increases systematically with p/, (see also the last 
column of Table [T] but note that the turnover masses there 
refer to the model GCMFs that we develop below). The sta- 
tistical significance of this is very high, and qualitatively it 
is the behavior expected if Mto owes its existence to cluster 
disruption at a rate that increases with p/,, as is the case with 
relaxation-driven evaporation. 

The right-hand panels of Figure |2] confirm once again that 
the GCMF peak mass is a very weak function of Galactocen- 
tric position. In fact, the observed distributions in the two rgc 



bins inside ~ 10 kpc are statistically indistinguishable in their 
entirety, and the main difference at larger rgc > 10 kpc is a 
slightly higher proportion of low-mass clusters rather than a 
large change in Mjo. All of this is consistent with the pri- 
mary dependence of the GCMF being that on ph, since Fig- 
ure[T]shows that the GC density distribution is not sensitive to 
Galactocentric position for rgc < 10-20 kpc but has a substan- 
tial low-density tail at larger radii (with a broader associated 
GCMF, as seen in the upper-left panel of Figure |2]). 

2.2. Simple Models 

We now assess more quantitatively whether these results 
are consistent with evaporation-dominated evolution of the 
GCMF from an initial distribution like that observed for 
young clusters in the local universe. We model the time- 
evolution of the distribution of M versus pi, in Figure[T]but do 
not attempt this for the distribution of ph over rgc — the details 
of which likely depend on a complicated interplay between 
the tidal field of the Galaxy, the present and past orbital pa- 
rameters of clusters, and the structural nonhomology of GCs. 
To compare our models to the cuiTent GCMF as a function of 
rgc, we simply calculate them using the observed ph of indi- 
vidual clusters in different ranges of Galactocentric radius. 

We assume that the initial GCMF was independent of clus- 
ter density, and that all globulars surviving to the present day 
have been losing mass for the past Hubble time at constant 
rates. We use the current half-mass density of each cluster 

1/2 

to estimate /iev « p^ . As we discussed earlier, an approxi- 
mately time-independent /igv is indicated by most calculations 
of two-body relaxation in tidally limited GCs. We give a more 
detailed, a posteriori justification in ^ for using p/,, rather 
than other plausible measures of cluster density, to estimate 

Pev- 

Consider first a group of coeval GCs with an initial mass 
function dN /d log Mq and a single, time-independent mass- 
loss rate /igv. The mass of every cluster decreases linearly as 
M(t) = M() - /ievf , and at any later time each has lost the same 
amount A = M^-Mit) = p^yt. FZOl show rigorously that in 
this case, the evolved and initial GCMFs are related by 

dN _ M dN _ M dN 

dlogM ~ Mo ^ d log Mo ~ (M+A) dlog(M+A) ' 

This is the basis for the claim that the mass function scales 
generically as dN /d log M oc M"^' (a /3 = power law) at low 
enough M(f ) < A — that is, for the surviving remnants of clus- 
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Fig. 2. — GCMF as a fun ction of half-m ass density, p/, = 3M/87rr,^ {left panels), a nd as a fun ction of Galactocentric radius, r^c (ri^ht panels), for 146 Milky 
Way GCs in the catalogue of lHarrisI 419961) . The dashed curve in all cases is an evolved lSchechteil function for the entire GC system jjordan et aU200'^ : equation 
(3) with 13 = 2, Mc = 10^ Mq, and A = 2.3 X lO' Mq for all clusters (from equation (4) and a median p), = 246 Mq pc"'), giving apeak atMxo = 1.6 X 10^ Mq. 
SoUd curves are the GCMFs predicted by equation (3) with fi = 2 and Mc = 10* Mq but individual A given by the observed pi, of each cluster (equation (4)) in 
the different subsamples. 



ters with Mq ~ A — just so long as the initial distribution was 

not a delta f unction . 

We follow FZOl' (see also lJordan et al.ll200 7^ in adopting a 
ISchech tei ( 1976) function for the initial GCMF: 

dN/d log Mo cx Mq"'^ exp (-Mq/M,) . (2) 

With /3 ~ 2, this distribution describes the power-law mass 
functions of young massive clusters in systems like the Anten- 
nae galaxies (e.g., Zhang & Fall 1999). An exponential cut- 
off at Mc ^ IQ Mq is generally consistent with such data, 
even if not always demanded by them; here we require it 
mainly to mat ch the curvature observed at high masses i n old 
GCMFs (e.g..l Burkert & Smithll2000t [Jordan et al.ll2007l) . 

Combining equations ([TJ and ^ gives the probability den- 
sity that a single GC with known evaporation rate and age has 
an instantaneous mass M. The time-dependent GCMF of a 
system of N GCs with a range of //gv (or ages, or both) is 
then just the sum of all such individual probability densities: 



dN 



d log M 



AiM 



[M+A,t 



exp 



M+A; 



Mc 



(3) 



Here the total mass losses A,- = {pL^^t)! may differ from clus- 
ter to cluster (f, being the age of a single GC) but both (3 and 
Mc are assumed to be constants, independent of ph in particu- 
lar.^ Given each A,, the normalizations A, in equation (O are 
defined so that the integral over d log M of each term in the 
s ummation is u nity. 

[Jordan et al.l (l2007l) have introduced a specialization of 
equation (O in which all clusters have the same A. They refer 

' Note that Mc appears to take on different values in the GCMFs of other 
galax ies, varying systematically with the total luminosity Lg;,] ( Jordan et al. 
I2007D . The reasons for this are unclear, as is the origin of this mass scale in 
the first place. 



to this as an evolved lSchechtej function and describe its prop- 
erties in detail (including giving a formula for the turnover 
mass Mto as a function of A and Mc) for the case (3 = 2. Here 
we note only that, at very young cluster ages or for slow mass- 
loss rates, such that A ^ Mc and only the low-mass, power- 
law part of the init ial GCMF is significantly eroded, any one 
evolved fSchechterl function has a peak at Mjo ~ A/(/3- 1) 
(for /3 > 1). As A increases relative to M^ the turnover at 
first increases proportionately and the width of the distribu- 
tion decreases (since the high-mass end at M > Mto is largely 
unchanged). For large A ;2> Mc, however, the peak is bounded 
above by Mto ^ Mc and the width approaches a lower limit. 
Thus, the dependence of Mto on A is weaker than linear when 
Mc is finite in the initial GCMF of equation (|2]l. Any peak in 
the full equation (O for a system of GCs with individual A 
values is an average of M different turnovers and must be cal- 
culated numerically. 

In their modeling of the Milky Way GC system, iFZOll ef- 
fectively compute mass functions of the type (O — ^based on 
the same initial conditions and dynamical evolution — with a 
distribution of A values determined by the orbital parame- 
ters of clusters in an idealized, spherical and static logarith- 
mic Galaxy potential (used both to fix ^ev in terms of clus- 
ter tidal densities and to estimate additiona l mass loss due 
to gravitational shocks). Jordan et af] ([2007) fit GCMF data 
in the Milky Way and scores of Virgo Cluster galaxies with 
their version of equation (O in which all GCs have the same 
A. They thus estimate the dynamical mass loss from typi- 

The increase of Mjq and the decrease of the full width of dN/d log M 
for increasing A eventually saturate when the mass los s per GC is so high 
that it affects clusters in the exponential part of the initial lSchechteil -function 
GCMF. This is because dN/d logM oc M+' exp(-M/Mc) is a self-similar 
solution to equation (TJ. 
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cal clusters in these systems. Here, we construct models for 
the MiUcy Way GCMF using A values given directly by the 
observed half-mass densities of individual GCs. 

We adopt /3 = 2 for the initial low-mass power-law index 
in equation (|2]l, whic h carries over into equation (O for the 
evolved dN/d log M. I Jordan et all t2007) have fitted the full 
Galactic GCMF with an evolved lSchechtei] function assuming 
[3 = 2 and a single A = A for all surviving globulars. They 
find ~ 10'' Mq and A = 2.3 x 10^ Mq. We use this value 
of Mc in equation (O and we associate A with the mass loss 
from clusters at the median half-mass density of the entire GC 
system, which is p), = 246 Mopc""* from the data in Figure [T] 

1 /2 

Since we are assuming that A = /ievf oc t for coeval GCs, 
we therefore stipulate 

A = 1 .45 X 1 0'* Mo {ph IMq pc-3 ) ' (4) 

for globulars with arbitrary p/,. Assuming a typical GC age of 
f = 13 Gyr, this corresponds to a mass-loss rate of 

/^ev ^ 1100 Mo Gyr-' (p,,/Mo pc'^) ''^ . (5) 

In ^ we discuss the cluster lifetimes implied by this value 
of /iev We emphasize here that the scaling of /igv and A 

1/2 

with Pfl follows rather generically from our hypothesis of 
evaporation-dominated cluster evolution, while the numerical 
coefficients in equations (|4|i and (|5]l are specific to the assump- 
tion of /3 = 2 for the power-law index at low masses in the 
initial GCMF. 

The dash ed curve s hown in every panel of Figure |2] is 
the evolved ISchechterj function fi tted to the entire GCMF of 
the Milky Way by iJordan et all (l2007l) . This has a peak at 
Mto ~ 1.6 X 10^ Mo (magnifiide My ~ -7.4 for a typical V- 
band mass-to-light ratio of 1.5 in solar units) and gives a very 
good description of the observed dN/d logM in the middle 
density bin, 79 < ph < 530 Mo pc"^, and in the two inner ra- 
dius bins, r„c < 9.4 kpc. This is expected, since the median 
half-mass density in each of these cluster subsamples is very 
close to the system- wide median p;, = 246 Mo pc"-' (see Ta- 
ble[T]i. Even in the outermost r^^ bin, a Kolmogorov-Smirnov 
(KS) test only marginally rejects the dashed-line model (at 
the ~95% level), because this subsample still includes many 
GCs at or near the global median pi, (see Figure [T]!. By con- 
trast, the average GCMF is strongly rejected as a model for 
the lowest- and highest-density GCs on the left-hand side of 
Figure 12] the KS probabilities that these data are drawn from 
the dashed distribution are < lO"'* in both cases. This is also 
expected since, by construction, these bins only contain clus- 
ters with densities well away from the median of the full GC 
system, for which the total mass lost by evaporation should be 
significantly different from the typical A = A(/5/,). 

The solid curves in Figure |2l which are different in ev- 
ery panel, are the superpositions of many different evolved 
ISchechten functions, as in equation (|3]l, with distinct A values 
given by equation (|4]i using the observed p/, of each cluster in 
the corresponding subsample. These models provide excel- 
lent matches to the observed dN/d log M in every pi, and rc,c 
bin, with x' < 1.3 per degree of freedom in all cases. This is 
the main result of this paper. 

The last column of Table [T] gives the mass Mjo at which 
each of the solid model GCMFs in Figure |2]peaks. We note 
that these turnovers increase roughly as Mto ^ P/f for 
our specific binnings in pi, and rgc, somewhat shallower than 

1/2 

the Pil scaling of the cluster mass-loss rate that defines the 



models. This is partly because of the averaging over indi- 
vidual turn overs implied by the summation of many evolved 
ISchechterl functions in each GC bin, and partly because — as 
we discusse d just after equation (O — the turnover mass of any 
one evolved lSchechted function cannot increase indefinitely in 

1 /^ 

direct proportion to A cx p^/', but has a strict upper limit of 
Mto <Mc. 

Our models are naturally consistent with the fact that the 
GCMF is narrower for clusters with higher densities. This 
is obvious in the left-hand panels of Figure |2j in the dis- 
cussion immediately after equation (O, we described how it 

1 /2 

follows from t he inc rease of Mto with A oc p^' for a sin- 
gle evolved ISchechte r function. In addition, the superposi- 
tion of many such functions with separate, density-dependent 
turnovers and widths results in wider GCMFs for cluster sub- 
samples spanning larger ranges of p/,. This accounts in partic- 
ular for the breadth of the mass function at rgc > 9.4 kpc. The 
globulars at these radii have 0.034 <Ph_<AAx 10^ Mo pc"^, 
corresponding to individual evolved lSchechter l functions with 
turnovers at 2.7 x 10^ < Mto < 4.0 x 10^ Mo. The compos- 
ite GCMF in the lower-right panel of Figure |2] is therefore 
extremely broad and shows a very flat peak, such that an over- 
all Mto cannot be establish ed precisely from the data alone. 
This explains the findings of iKavelaars & HanesI (Il997i) . who 
pointed out that the GCMF of the outermost third of the Milky 
Way cluster system has a turnover that is statistically consis- 
tent with the full-G alaxy average, but a larger dispersion (see 
also lGnedinlll997h . 

Finally, if the GCMF evolved dynamically from initial con- 
ditions similar to those we have adopted, then the data and 
models in the left-hand panels of Figure |2] argue against 
the notion that external gravitational shocks, rather than in- 
ternal two-body relaxation, were primarily responsible for 
shaping the present-day GCMF. This is because the mass- 
loss rate caused by shocks alone, -dM/dt = psh M / ph, 
differs significantly from that caused by evaporation alone, 

1 /2 

-dM/dt = pev oc Pil . The direct dependence of psh on M en- 
sures that shocks become progressively less important com- 
pared to evaporation as clusters lose mass (at a given p/,), and 
consequently shocks are not likely to have had much effect on 
the observed GCMF for M < Mtq. Furthermore, the inverse 
dependence of psh on p/, is contrary to the direct dependence 
of Mto on p/, shown in Figure |2] The different roles played 
by shocks and evaporati on in s haping the observed GCMF are 
discussed more fully bv lFZOll We note here that gravitational 
shocks may have been important in destroying very massive 
or very low-density clusters early in the history of our Galaxy. 

2.3. Other Cluster Properties 

If the current shape of the GCMF is fundamentally the re- 
sult of long-term cluster disruption according to a mass-loss 

1 /2 

rule like pev oc Pi/ , then it should be possible to reproduce 
the distribution as a function of any other cluster attribute by 
using the observed p/, of individual GCs in equations (O and 
(|4|i to build model dN/d log M for subsamples of the Galactic 
cluster system defined by that attribute — as we did for the r^^ 
binning of ^2.21 Here we explore one example in which dif- 
ferences in the GCMFs of two groups of globulars can be seen 
in this way to follow from differences in their p/, distributions. 

Smith & Burkert (2002) have shown that the mass function 
of Galactic globulars with King ( 1966) model concentrations 
c < 0.99 has a less massive peak than that for c > 0.99. [Here 
c = log(r,/r()), where r, is the fitted tidal radius and ro a core 
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scale.] They further find that a power-law fit to the low- 
c GCMF just below its peak returns dN /d logM oc M"^"^ — 
shallower than the M"""' expected generically for a mass-loss 
rate that is constant in time — but they confirm that the latter 
slope applies for the GCMF at c > 0.99. They discuss various 
options to explain these results, including a suggestion that, if 
the mass functions of both low- and high-concentration clus- 
ters evolved slowly from the same, young-cluster-like initial 
distribution, then the mass-loss law for low-c GCs may have 
differed from that for high-c clusters. However, they give no 
physical explanation for such a difference, and we can show 
now that none is required. 

The upper panel of Figure |3] plots concentration against 
half-mass density for the same 146 GCs from Figure [T] the 
filled circles distinguish 24 clusters with c < 0.99. There 
is a correlation of sorts between c and ph, which either de- 
rives from o r causes the better-known correlation between c 
and M (e.g., 'Djorgo vski & Meylanlll994t lMcLaughUnll2000l) . 
The important point here is that the pf, distribution is off- 
set to lower values and has a higher dispersion at c < 0.99. 
Following the discussion in ^2.21 we therefore expect the 
low-concentration GCMF to have a smaller Mto, a flatter 
shape around the peak, and a larger full width than the high- 
concentration GCMF. 

The lower panel of Figure|3]shows the GCMFs for c < 0.99 
(filled circles) and c > 0.99 (open circles). The curves are 
again given by equation (O with (3 = 2, = 10^ Mq, and in- 
dividual A calculated from the observed cluster p/, through 
equation (|4]l. These models peak at Mto — 4.3 x 10"^ M© 
for the c < 0.99 subsample but at Mjo ~ 1.8 x lO'' Mq for 
c > 0.99, entirely as a result of the different pi, involved. 
The larger width of dN / d log M and its shallower slope at 
any M < 10^ Mq for the low-concentration GCs are also 
clear, in the model curves as well as the data. It is further 
evident that there are no low-c Galactic globulars observed 
with M > 2 X 10^ Mq, above the n ominal turnover of the full 
GCMF (as lSmifli & Burkertll2002l noted). But this is not sur- 
prising, given that there are so few low-concentration clusters 
in total and they are expected to be dominated by low-mass 
objects because of their generally low densities. Thus, the 
solid curve in Figure |3] predicts perhaps ~3 high-mass clus- 
ters with c < 0.99, where none is found. 

The apparent variation of the Milky Way GCMF with in- 
ternal concentration is therefore consistent with the same 
density-based model for evaporation-dominated dynamical 
evolution that we compared to dN /d log M as a function of 
Pi, and rgc in ^2.21 To show this, we have made use of the 
densities pi, exactly as observed within the two concentration 
bins indicated in Figure [5]— just as we also took pi, directly 
from the data for GCs in different ranges of rgc to construct 
models for comparison with the observed dN /d log M in the 
right-hand panels of Figure |2] Of course, this is not the same 
as explaining the distribution of p/, versus rgc or c. Doing so 
would certainly be of interest in its own right, but it is beyond 
the scope of our work here. 

3. DISCUSSION 

In this section, we first show that the mass-loss rate in equa- 
tion (|5]l above implies cluster lifetimes that compare favorably 
with those expected from relaxation-driven evaporation. Then 

1/2 

we discuss why it is reasonable to approximate /igv Pi, in 
the first place. Finally, we address the issue of possible con- 
flict, in some other models for evaporation-dominated GCMF 
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Fig. 3. — Top: Concentration parameter as a function of half-mass density 
for 146 Galactic GCs. The line of points at c = 2.5 comes from the practice of 
assigning this value to core-collapsed clusters in the Harris (1996) catalogue 
and its sources. Bottom: GCMF data and models (eqs. (5] and (4)) for 24 
clusters with c < 0.99 (filled circles and solid curve) and 122 clusters with 
c > 0.99 (open circles and dashed curve). 

evolution, between the near-constancy of Mjq as a function 
of rgc and the observed kinematics of GC systems. 

3.1. Cluster Lifetimes 

The disruption time of a GC with mass M and a steady 
mass-loss rate /iev is just fjis = Mj p^^,. It is convenient, 
for purposes of comparison with evaporation times in the 
literature, to normalize fdis to the relaxation time of a 
cluster at its half-mass radius. In general, this is f,h = 
0.138M'/V^/^/ [G'/^m*ln(7M/mj], where is the mean 
stellar mass. For clusters of stars with a single mass , 
c± O.IMq, and j .4 are appropriate dSpitzeii 119871 : 
iBinney & Tremain3 Il987l equation [8-72]), in which case 
equation Q for /igv from our GCMF modeling implies 



fdis 
tih 



M 

Mev^rh 



10 



In {0.51 M/Mq 
In (0.57 X 105) 



(6) 



Clusters with realistic stellar mass spectra will have slightly 
different values of m^, and a smaller 7 in t he calculation of 
the relaxation time ( Giersz & Hegg ie 1996), which changes 
the numerical value of fdis/frh somewhat but does not alter any 
scalings. 
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We obtained the normalization of /iev oc p^' in ^2.21 by 
fitting to observed GCMFs constructed by applying a spe- 
cific mass-to-light ratio Ty to every cluster, with models as- 
suming a specific form for the initial dN/d logMo. Thus, 

the result in equation (|6j depends both on the median Ty 
and on the power-law index f3 at low masses in the original 
ISchechterl -function GCMF. The net scaling, for either single- 
or multiple-mass clusters, is 

fd,sAi-hCxT-'/^/3-l)-i . (7) 

To see the dependence of this dimensionless lifetime on 
Ty, note that we require /iev oc A cx Ty to fit the mass losses 
of clusters with a given distribution of luminosities (the di- 

1 /2 

rect observables), whereas M/fih is proportional to p,/ oc 

Tl/\L/riy/\ Therefore, fd.s/frh « (M/f,h)/Mev « Ty''^ The 
mass-to-light ratios adopted in this paper, with a median value 

Tv ~ 1.5 Mq Lq, are tied directly to dynamical determina- 
tions ( ET i. 

To understand the dependence on (3 in equation d?), recall 
first that the coefficients in our expressions for A and /igv 
as functions of pi, (eqs. H] and ||5|) followed from choos- 
ing P = 2 for the power-law exponent at low masses in the 
initial GCMF (equation ||2l). As we mention ed just after 
equation (O, the turnover mass of an evolved Schechtgi func- 
tion with any /3 > 1 is Mjq ~ A/(/3- 1) in the limit of low 

1 /2 

A OC , and Mjo — > Mc for very high A. In this sense, 
the strongest observational constraints on the normalizations 
of A and /iev come from the low-density clusters. All other 
things being equal, their GCMF can be reproduced with 13 2 
if A and /iev are multiplied by (/3- 1) at fixed p/,. Therefore, 
fdis oc 1 // Hev OC 1/(0- 1). Obse rvations of young massive clus- 
ters (e.g.. lZhang & Fallll 19991) indicate that /3 is near 2; but if 
it were slightly shallower, then the cluster lifetimes we infer 
from the old GCMF would increase accordingly. Even a rela- 
tively minor change to /3 = 1 .5 would double fdis/frh from « 10 
to «20. 

In the model of iHenonI (Il961h for single-mass clusters 
evolving self-similarly (fixed ratio pt/ph of mean densities in- 
side the tidal and half-mass radii) in a steady tidal field (p, 
constant in time), a cluster loses 4.5% of its remaining mass 
every half-mass relaxation time. The time to complete disrup- 
tion is therefore fdis/frh = 1 /0.045 ~ 22. For non-homologous 
clusters in a steady tidal field, fdis/frh is a function of central 
concentration and can differ from the Henon value by factors 
of abo ut two. From one-dimen sional Fokker-Planck ca lcula- 
tions. iGnedin & Ostriked (11997 ) find fdis/fih ^ 10-40 for|King| 
d 19661) model clusters with c values similar to those found in 
real GCs and with gravitational shocks suppressed (see their 
Figure 6). Thus, even though the evaporation time in equation 
Q may be slightly shorter than is typically found in theoreti- 
cal calculations, it is certainly within the range of such calcu- 
lations. Moreover, the assump t ions o f a st eady tidal field and 
a single stellar mass in 'H enonl(ll961h and lGnedin & Ostrikerl 
(jT997) are important. Part of the difference between the typ- 
ical lifetimes in these particular theoretical treatments and 
our estimate of fdis/fih from the GCMF is that the former do 
not include gravitational shocks, which may have accelerated 
somewhat the evolution of real clusters (although we stress 
again that shocks do not appear in general to have dominated 
the evolution of extant Galactic GCs and are not expected 
to affect the basic time-independence of the net mass-loss 
rate; see IVesperini & Heggi3ll997l iGnedin. Lee. & Ostriked 



ll999llFZ0lL and lPrieto & Gnedinll2006l) . A spectrum of stel- 
lar masses in the clusters may also have contributed to an in- 
crease i n evapor ation rate over th e single-mass values (e.g., 
iJohnsto ne 1993; Lee & Goodman|[T9 95). 

Estimates of evaporation times from other numerical 
methods and for models of multimass clusters can be rather 
sensitive to the detailed computational techniques and input 
assumptions and approximations, and differences at roughly 
the factor-of-two level in fdis/fit ] between different analyse s 
are not uncom mon; see, e.g ., Vesperini & Heggie ( 1997^ 
iTakah ashi & Po rtegies Zwarll (1998 200(3) ' Baumgardl 
(2001), Joshi. Nave. & Rasiol Hooih . lOierszl md 
Baum gardt & Makinol (l2003h . Thus, although the lifetimes 
in these studie s ten d to be broadly comparabl e to those in 
iHenonI ( 1196 lb and IGnedin & Ostrikei^ (1 19971) . noticeably 
shorter values do occur in some models. In any case, we 
are encouraged by consistency to within factors of two or 
three between estimates of fdis or /iev by such vastly different 
methods — one purely observational, based on the mass 
functions of cluster systems; the other purely theoretical, 
based on idealized models for the evolution of individual 
clusters — particularly since each method involves several 
uncertain inputs and parameters. 

3.2. Approximating p^v p]/^ 
3.2. 1 . Half-mass versus Tidal Density 

The dimensionless disruption time in equation ^ is inde- 
pendent of any cluster property other than the Coulomb log- 
arithm because we have used GC half-mass densities to esti- 

1/2 1/2 

mate fdis = M//Xev oc M/p^' , while fi-h also scales as M/ pjl . 
However, as we mentioned above, t he Fokker-Planck calcu- 
lations of Gnedin & Ostri keJ (11997 ) in particular show that 
fdis/fi-h is ac tually a function of central concentration, c, for 
lKingl(ll966l) model clusters in steady tidal fields. The constant 

1/2 

of proportionality in /iev should therefore also depend 

on c, a detail that we have neglected to this point. We show 
now that this has not biased any of our analysis or affected our 
conclusions. 

The dotted curve in Fig ure HI i llustrates the dependence of 
fdis/frh on c for singlc-mass King models, as given by equation 
(30) of Gnedin & Ostriker ( 1997). The solid curve is propor- 
tional to (ph/ PtY^'^ = (rf /2riy^^, which we have calculated 
as a function of c for these models and multiplied by a con- 
stant to compare directly with fdis/frh- Evidently, there is an 
approximate equality fdis/frh ~ 2.15(p/,/ ptY^', which holds to 
within < 15% over the range of concentrations shown in Fig- 
ure |4] (note that all but 6 Galactic GCs have 0.7 < c < 2.5, 
corresponding to central potentials 3 < Wo < 1 1). Thus, if the 

evaporation time is written as fdis oc t^hiph/ Pt)^^^ o^Mj pj'^, 
then the constant of proportionality in the mass-loss rate /iev oc 

M/fdis oc pj'^ should be nearly independent of c. In fact. lKingI 
(11966) originally concluded, from quite basic arguments, that 
the evaporation rate of a cluster with a lowered-Maxwellian 

1 11 

velocity distribution would take the form /iev oc p, with 
only a weak dependence on c. An essentially concentration- 

1/2 

independent scaling of /iev with p, is also found in A^- 
body si mulations of tidally limited, multimass clusters (e.g., 
IVesper ini & Heggie 1997) and so is not an artifact of any as- 
sumptions specific to the calculations of either lKingI (1 1 9661) or 
.Gnedin & Ostriker (1997.) . 

This suggests that it might have been more natural to spec- 
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c = log(ryrJ 



Fig. 4. — Dependence of ?dis/?rh (dotted line; from IGnedin & Ostrikeil 
119971) and (p^/pr)'^^ (solid line; after scaling by a factor of 2.15) on cen- 
tral concentration for single-mass King-model clusters. Over the range of c 
shown, which includes nearly all Galactic globulars, the approximate propor- 
tionaUty fdis/frh °^ iph/Pi)^^^ holds to within better than 15%. Thus, to this 
level of accuracy the evaporation time tus is roughly the same multiple of 
1 /2 

M/p,' for clusters with any internal density profile. 

1/2 

ify cluster evaporation rates proportional to p, rather than 

p]/'^ when developing our GCMF models in ^ For any clus- 
ter in a steady tidal field, with a constant p,, such a choice 
would also have been automatically consistent with an ap- 
proximately time-independent /igy and the corresponding lin- 
ear M(f) dependence that we have adopted throughout this 
paper As we discussed at the beginning of ^ our deci- 
sion to work with pi, rather than p, was motivated by the 
fact that the half-mass density is much better defined in prin- 
ciple and more accurately observed in practice. Neverthe- 
less, re-writing oc pj^^ as ^ev oc {pi/ phY^^ x p]/^ makes 
it clear that the validity of our models, with a fixed coefficient 

1/2 

in yUev oc Pf^ , depends on the extent to which variations in 
(pt/phy^'^ can safely be ignored. 

Figure H] shows that the full range of possible values for 
(ph/ ptY^^ in lKing '-model clusters with c > 0.7 is only a factor 
of ~4 between minimum and maximum. Therefore, using 
a single, intermediate value of this density ratio to describe 
all GCs (or a single GC evolving in time through a series of 
quasi-static King models) — which we have effectively done 
by using a GCMF fit to normalize A and p^^y in equations (|4|i 
and (|5]l — should never be in error by more than a factor of 2 or 
so. This is a relatively small inaccuracy, given that measured 
GC densities range over four to five orders of magnitude. 

1/2 

To confirm more directly that our models with /iev p,/ 
are good approximations to GCMF evolution under a mass- 
loss law /iev oc pj'^ , we have repeated the analysis of ^in full 
but using t he GC tidal d ensities pt (derived from the values of 
r, listed bv Harri sill 9961) in place of p/, throughout. All of our 
main results persist. 

For example, the two panels of Figure |5] which are analo- 
gous to the left- and rightmost panels of Figure[T]above, show 
that (1) the GC mass distribution has a clear dependence on 



Pt, with a lower envelope that is well matched by a line of 

1 /2 

constant evaporation time, M (x p,' (the dashed line in the 
plot); and (2) although the scatter in the distribution of p, 
over Galactocentric radius is smaller than the scatter in p/, ver- 
sus Tgc, it is still significant. Because the M-rgc distribution 
can now be viewed as the convolution of the M-p, distribu- 
tion with the Pt-fgc distribution, the scatter in p, versus r^c is 
again critical in explaining the weak or null dependence of the 
GCMF on Galactocentric radius. (The M-r„c distribution is, 
of course, unchanged from that shown in the middle panel of 
Figure!!])" 

Figure|6]shows the Milky Way GCMF for globulars in three 
equally populated bins of tidal density (defined as indicated in 
the left-hand panels of the plot) and in the same three bins of 
Galactocentric radius that we used in ^2. 21 above. Our mod- 
els for these distributions are based as before on equation (O 
with (5 = 2, but now the total mass lost from any GC is esti- 
mated from its tidal density rather than its half-mass density. 
Specifically, we take 

A = 2.1 X 10^ Mq (p,/M0pc-^)'''^ . (8) 

The numerical coefficient in equation (O is such that it gives a 
A identical to that in equation ^ for a GC with p/,/ p, = 210, 
whic h is the m edian value of this density ratio for the 146 GCs 
in the lHarrisI ( IT996.) catalogue. 

As in Figure |2l the dashed curve in every panel of Figure 
|6]is the same, representing a fit to the average dN/d log M of 
the entire Galactic GC system. Thus, it is immediately clear 
that the peak mass of the GCMF increases significantly and 
systematically with increasing p,, just as it does with increas- 
ing Pi,. Meanwhile, the solid curves are subsample-specific 
model GCMFs, obtained by using the observed tidal density 
of each cluster in any p, or rgc bin to specify i ndividual A val- 
ues via equation ([8]) for each of the evolved Schechter func- 
tions in the summation of equation As expected, there is 
no appreciable difference, in terms of the fits to any of the ob- 
served GCMFs, between these models based on evaporation 

rates oc p^^ and our original models with p,ev oc p'/^. 

3.2.2. Retarded Evaporation 

Another po tential concern comes from recent arguments 
(see especiallv lBaumgardtl2001l:lBaumgardt & Makinol2003h 
that the total evaporation time of a tidally limited cluster 
is not simply a multiple of an internal two-body relaxation 
time, frix oc (Mr-*)'/^, but depends on both fi-ix and the crossing 
time fcr oc (M/r^)~^^^ through the combination fjis oc t^^J^^^, 
with jc < 1. The mass-loss rate p,ev ocM/fdis then scales as 
^3/2-Ay-3/2^ which for x=/ I differs from the rates p,ev oc p]/^ 

and /iev oc p^^ that we have so far adopted. However, our 
GCMF models are still meaningful, because postulating fdis oc 
^rix^cr"^ implies a dependence of /iev on a measure of cluster 

1/2 

density that is, once again, well approximated by p^^ for 
Galactic GCs. Before showing this, we briefly discuss the 
reasons and the evidence for a possible dependence of fjis on 
both fi-ix and t^. 

' ' As was also the case with our earlier plots involving p;, in Figure \T\ 
the scatter and structure in both panels of Figure \5\ are real, since the rms 
scatter of log r, about the best-fit lines to either of log M or log rgc is 0.3-0.35 
while the rms errorbars based on formal fitting uncertai nties are in the range 
(5(log r,) ~ 0.05-0.15 for a variety of models tMcLaughlin & van der Mm^ 
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p, [M^ pc-3] r^, [kpc] 

Fig. 5. — Scatter pl ots of mass M versus mean density inside the tidal radius (p, = 3M/47rr,') and ofp, versus Galactocentric radius Tgc, for 146 Galactic GCs 

from the IHarri3 419961) catalogue. These plots are analogous to the left- and rightmost panels of Figure[T] The dashed line in the left-hand plot traces the relation 

1/2 1/2 
M oc p/ , which defines a locus of constant evaporation time for p^v on p/ . 
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Fig. 6. — Observed GCMF (points, with Poisson errorbars) and models (curves) as a function of mean cluster density inside the tidal radius, p, = 3M/4iTrf 
(left-hand panels), and as a function of Galactocentric radius, rgc (right-hand panels). The dashed curve in every panel is an evolved Schechter function 
representing the entire GC system: equation (3) with /3 = 2, Mc = 10* Mq, and a single A, common to all clusters, evaluated from equation (8) using the median 
Pi of all 146 Galactic GCs. Solid curves are subsample-specific models using equation (3) with /3 = 2 and Mc = 10* Mq but a different A value for every cluster 
(obtained from equation (8] using individual observational estimates of p,) in any p, or fgc bin. 
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If stars are assumed to escape a cluster as soon as they 
have attained energies above some critical value as a result 
of two-body relaxation, then fdis oc frix is expected (a nd con- 
firmed by ^-body simulations; e.g. jBaumgardtilToOlh . How- 
ever, more complicated behavior may arise when escape not 
only depends on stars satisfying such an energy criterion, but 
also requires them to cross a spatial boundary. Then, al- 
though the stars are still scattered to near- and above-escape 
energies on the timescale f^ix, they require some additional 
time to actually leave the cluster This escape timescale is 
related fundamentally to (but also depends on details of 
the stellar orbits, the external tidal field, and the shape of 
the zero-energy surface). The longer this extra time, the 
higher is the probability that further encounters with bound 
cluster stars may scatter any potential escapers back down 
to sub-escape energies. The net result is a slow-down ("re- 
tardation") of the overall evaporation rate ([Cha ndrasekhar 
1942"; 'Kins"1959^, 'Takahashi & Portegies Zward[1998, ,2000; 
Fukushige & Heggie 2000; Baumgardt 200 iT and a length- 
ening of the cluster lifetime fdis, by a factor that can be ex- 
pected to increase with the ratio fa /frix- If this factor scales as 

(fcr/frlx)'"-* for some X < 1, then fdis oc frlx(fci/fi-lx)'"' = t-^Ja"- 

While such a retardation of evaporation can be expected 
to occur at some level in all clusters, there are physical sub- 
tleties in the effect that are probably not captured adequately 
by a simple re-parametrization of lifetimes as fdis oc t^^j^'". 
In particular, it is unlikely that this expression can hold for 
clusters of all masses with a single value of ;ic < 1. Since 
fcr/fiix <x M"', very massive clusters have fcr ^ frix, and stars 
scattered to greater than escape energies by relaxation cross 
the tidal boundary effectively instantaneously — implying that 
the standard fdis oc fiix, orx ^ 1, applies in the high-mass limit. 
Indeed, if this were not the case, and a fixed x < 1 held for all 
M, then an unphysical ? d is < ^r ix would obtain at high enough 
masses; see Baumg ardt! (1200 ll) for further discussion. Unfor- 
tunately, "very massive" is not well quantified in this context, 
and it is not yet clear if a single value of x is accurate for the 
entire GC mass regime. So far, it has been checked directly 
only for initial cluster masses below the current peak of the 
GCMF. 

It is also worth noting that the analysis and simulations 
aimed at this problem to date have dealt with clusters on cir- 
cular or moderately eccentric orbits in galactic potentials that 
are static and spherical. This means that any tidal perturba- 
tions felt by stars within the clusters are relatively weak and/or 
slow compared to their own orbital periods, leading to nearly 
adiabatic or at least non-impulsive responses. In more realis- 
tic situations, the galactic potential would be time-dependent 
and non-spherical and there might be additional tidal pertur- 
bations, including disk and bulge shocks. These perturbations 
could in some cases accelerate the escape of weakly bound 
stars from the clusters and thus counteract the retardation ef- 
fect to some degree. Further study is therefore needed to de 
termine the regime of validity of the formula fdis oc f, 
its possible modification outside this regime. 

In the m eantime, Baumga rdt ((200lT) and 
iBaumgardt & Makinol (120031: hereafter IBM03I) have fit- 



rixC' and 



ted this formula to the lifetimes of a suite of A^-body clusters 
with initial masses Mn 5, 7 x Ify^Mp, and several different 
initial concentrations and orbital eccentricities. BM03' at 
first write fdis in terms of the relaxation and crossing times 
of clusters at their half-mass radii, so that fiix oc (Mr\yi^, 
fcr OC (M/r^)"'/^, and fdis oc M^'^l'^r]^^ (see their equation 



[5]). However, they immediately take a factor of {rtjrhY'l'^ 
out from the normalization of this scaling — in effect to obtain 

fdis oc M^~^I^PJ^ with a different constant of proportionality — 
and then use a simple definition of the tidal radius (their 
equation [1], = GMrp/2V^, which is appropriate for a 
circular orbit of radius in a logarithmic potential with 
circular speed Vc', see llnnanen. Harris. & WebbinS Il983h 
to obtain the total lifetime of a cluster as a function of its 
initial mass, perigalactic distance, and Vc (their equation 
[7]). A single exponent x ~ 0.75 and a single normalization 
in this function then suffice to predict to within 10% the 
lifetimes of the simulated clusters, regardless of their initial 
concentrations. By implication, if frix and f^r were fixed at r/, 
rather than r,, then fdis would have an additional concentration 
dependence, related t o the ratio {r,/rh)^l^ — very similar to 
what we discussed in ^3.2.1l for the case x=\. 

We now re-examine the Milky Way GCMF in terms of 
this prescription for retarded evaporation (bearing in mind the 
caveats mentioned above). To avoid any explicit dependences 
on concentration, we also focus on the tidal radius and write 
fdis oc M^~^l^r]^^ for general x < 1; but we do not substitute a 
potential- and orbit-specific formula for r, in terms of rp and 
galactic properties such as Vc- Instead, to keep the empha- 
sis entirely on cluster densities, we re-write the scaling of the 
lifetime in terms of the mean surface density inside the tidal 
radius, S, = M/irr}, and the corresponding volume density 



p, = 'iM/Airr'^. This leads to fdis oc MS, 
then implies 



3(l-x) -2(x-3/4) 
Ft 



^ev = -dM/dt OC M/fdis OC Sf<'"-'V?*'' ^^"^^ 



I /2 

Clearly, the standard ^ev oc p, , which we have already dis- 
cussed, is recovered for x=l; while for x = 0.75, we have the 

3/4 

e qually s traightforward fi^v oc E, . 

|BM03| find that, even with the retarded evaporation implied 
by X ~ 0.75, the masses of their simulated clusters still de- 
crease approximately linearly with time after stellar-evolution 
effects (which are only important for the first few 10** yr) are 
separated out; see especially their Figure 6, equation (12), and 
related discussion. Thus, if the GCMF initially rose towards 
low masses and has been eroded by slow, relaxation-driven 
cluster destruction, then in this modified description of evap- 
oration we might expect the current mass function to depend 
fundamentally on E, rather than pi, or p,. But because M(f) 

still decreases nearly linearly with f , only now with /igv oc T,'^^'^ 
for each cluster, the shape of the evolved GCMF and its de- 
pendence on E, should resemble our earlier results for /?/, and 

Pi- 

We have confirmed this expectation by repeating all of our 

analyses in ^again, now using oc E,^'''* to estimate cluster 
ma ss-loss r ates. As before, we calculate E, from the data in 
the iHarrisI (1 19961) catalogue, although we caution once more 
that the tidal radii, and thus the derived E,, are more uncertain 
than ri, and p/,. 

Figure |7] which should be compared to Figures [T] and |5] 
above, shows that the average Galactic GC mass increases 
systematically with E, ; that the lower envelope of the M-E, 

3/4 

distribution is described well by M oc E, (the dashed line in 
the left-hand panel of Figure which is a locus of constant 

3/4 

lifetime against evaporation for /igv oc E, ; and that the scat- 
ter in the distribution of cluster E, versus Galactocentric ra- 
dius (right-hand panel of the figure) is substantial, as required 



, which 



(9) 
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to account for the almost non-existent correlation between M 
and Tgc . 

The left-hand side of Figure [8] shows the mass functions 
of globulars in three bins of S,, as defined in each panel. 
The right-hand side of the figure shows dN /d log M in the 
same three intervals of rgc as in Figures |2] and |6] above. As 
in those earlier plots, the dashed curve in all panels of Figure 
[8] is a model GCMF with the same parameters in every case, 
representing the mass function of the entire Galactic GC sys- 
tem. Once again, compared to the average Mjo, the observed 
turnover mass is significantly lower for clusters in the lowest 
S, bin and higher for clusters in the highest E, bin, while the 
width of dN /d log M decreases noticeably as S, increases. 

The solid curves in Figure |8] are again differe nt in every 
panel. They are the sums of evaporation-evolved ISchechterl 
functions as in equation (O, with the usual (3 = 2 assumed 
but with total mass losses estimated individually for each GC 

in any E, or rgc bin according to A ex T^J'^ rather than A oc 

1/2 1/2 

Pfl or A OC p, . However, it turns out not to be necessary 

to change the normalization of A oc p]^^ in equation (|4|i to 
achieve good fits to the observed GCMF as a function of either 
E, or Tgc. Thus, in Figure|8]we have simply used 



A = 1.45 X IO'^Mq (E,/M0pc- 



^3/4 



(10) 



The fits of these models, based on fjis oc t'^xjl'^" with x ~ 
0.75, are indistinguishable from the fits of our original mod- 
els based on the standard fdis oc fiix, i.e., x = 1. (We have con- 
firmed that adopting individual A given by equation ifTOl also 
reproduces the GCMFs of low-and high-concentration GCs 
in Figure [3] as well as before.) It was somewhat unexpected 
that equation ( fTOl i and equation (|4|i should have the same nu- 
merical coefficient, but we note that this follows empirically 
from the fact that the measured ph and E, of Galactic GCs 
are consistent with the simple near-equality, pi,/MqPC~^ w 
(E,/M0pc-2)i-5 in the mean. This is illustrated in Figure |9] 
which also shows that there is significant scatter about the re- 
lation.'^ However, this scatter does not correlate with clus- 
ter mass or Galactocentric radius. From a pragmatic point of 

view, therefore, and T,^^^ are near enough to interchange- 
able for our purposes, and there is no practical difference be- 
tween GCMF models based on one or the other measure of 
GC density. 

One further check on this is to verify that the mass-loss rate 
associated with equation (fTOl i is roughly in keeping with that 
implied by the A^-body simulations pointing to jc = 0.75 in the 
first place. Thus, we compare the rate 

= A/(13 Gyr) ~ llOOMeGyr"' (E,/Mopc-2)3/4 (H) 

to a formula implicit in IBM03L Starting with their equa- 
tion (7) for the lifetime fjis as a function of initial cluster 
mass and perigalactic distance and circular speed in a loga- 
rithmic halo potential; using their x = 0.75 and their normal- 
ization of 1.91 X 10^ yr, multiplied as in their equation (9) 
by ( 1 + e) to allow for eccentric orbits with apo- and peri- 
galactic distances related by e = (rg- rp)/{ra + rp); insert- 
ing their equation (1) for r,; taking the mean mass of clus- 
ter stars to be = O.55M0, as they do; using 7 = 0.02 as 

Although it may be only a coincidence that the constant of proportion- 
ality in Pi, oc S/^ is so near unity, the basic scaling itself holds because 
combining the observed correlation between cluster mass and central con- 
centration liDiorgovski & Mevlan. 1994; McLauphlin. 2000.) with the intrinsic 
dependence of r, /ri, one in lKind models leads roughly to (r, /r^ ) oc M' . 



they do in the Coulomb logarithm, ln(7M()/m,); and defining 
E, = Mo/ttt^q (the subscript denoting initial values), we 
obtain 



^ev(BM03): 



0.7Mo 



560 ^ _, 



fdis 1+e 

"in (O.O36MO/M0 
In (0.036 X 105) 



3/4 



Mopc- 



3/4 



(12) 



This is appropriate for clusters that just fill their Roche lobes 
at perigalacticon, which is where E,.o is specified. The factor 
of 0.7 in the first equality accounts for mass loss due to stellar 
evolution in the BM03 simulations, which, as they discuss, 
can be treated as having occurred almost immediately and in 
full at the beginning of a cluster's life. 

Our GCMF-based /igv is a factor of « 2 faster than the A^- 
body value for clusters on circular orbits (with e = and in 
steady tidal fields) in the simulations; and our fi^y is still 
within a factor of about three of the A^-b ody rate for clusters on 
eccentric orbits with e = 0.5 in lBM03l (e ~ 0.5-0.6 is typical 
for tracers wit h an isotropic velo city distribution in a logarith- 
mic potential: Ivan den Bosch e t al . 199 9). This is very similar 
to the comparison of lifetimes in ^3.1l for our original models 

based on /iev oc . Moreover, our new estimate of /lev and 
that in BM03 are still subject to their own, separate uncertain- 
ties and reflect different idealizations and assumptions. For 
example, our rate still depends on the exact power-law expo- 
nent f3 at low masses in the initial GCMF, as discussed after 
equation (|7]); while the rate from BM03 still neglects grav- 
itational shocks from disk crossings and passages by a dis- 
crete galactic bulge, and may additionally be biased low for 
Mo > lO^M© if X > 0.75 at such masses. All of this — not 
to mention again the large uncertainties and possible system- 
atics in the estimates of tidal radii needed to calculate E, — 
makes the near agreement between equations ( fTTI ) and ( fT2] i 
more striking than any apparent discrepancy. 

In summary, although the relation p.. 
in time is rigorously correct only in rather specific circum 
stances, our GCMF models based on it in ^are good proxies, 
in all respects, for models based on other plausible characteri- 
zations of relaxation-driven cluster mass loss. This result will 
likely be important for future studies of the mass functions of 
extragalactic cluster systems, where it may well be necessary 
to adopt procedures based on pi, rather than p, or E, because 
of the difficulty or impossibility of estimating tidal radii. 

3.3. Mto versus rgc, and Velocity Anisotropy in GC Sy stents 

In this paper we have directly modeled dN /d log M as a 
function only of GC density and age, and used the observed p/, 
(or pf, or E,) of clusters in relatively narrow ranges of Galac- 
tocentric position to show that such models are consistent with 
the current near-constancy of the GCMF as a function of r^^- 
Most other models in the literature for evaporation-dominated 
GCMF evolution, in either the Milky Way or other galaxies, 
instead predict the distribution explicitly as a function of rgc at 
any time. They therefore need, in effect, to derive theoretical 
density-position relations for clusters in galaxies alongside 
their main GCMF calculations. This usually begins with the 
adoption of analytical potentials to describe the parent galax- 
ies of GCs. Taking these to be spherical and static for a Hub- 
ble time allows the use of standard tidal-Umitation formulae 
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Fig. 7. — S catter plots of mass M versus mean surface density inside tlie tidal radius (S, = M/irrf) and of E, versus Galactocentric radius Cgc, for 146 Galactic 

GCs from the lHarrisI 119961) catalogue. These plots are analogous to the left- and rightmost panels of Figure[T] and the two panels of Figure|5] The dashed line in 

^/4 3/4 
the left-hand plot traces the relation M ocT,, , which defines a locus of constant evaporation time for fi^y oc S, . 
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Fig. 8. — Observed GCMF (points, with Poisson errorbars) and models (curves) as a function of mean surface density inside the tidal radius, E, = M/nrf (left- 
hand panels), and as a function of Galactocentric radius, fgc (right-hand panels). The dashed curve in every panel is an evolved Schechter function representing 

the entire GC system: equation ^3) with (3 = 2, Mc = 10* Mq, and a single A, common to all clusters, evaluated from equation \10\ using the median S, of all 
146 Galactic GCs^ Solid curves are subsample-specific models using equation (3) with /3 = 2 and Mc = 10* Mq but a different A value for every cluster (obtained 
from equation IIQI using individual observational estimates of E,) in any or rgc bin. 
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Fig. 9. — Half-mass density, pi, = SM/Sirr^^, against mean surface den- 
sity inside tlie tidal radius, = M/nrf, for 146 clusters with data in lHarrisI 
(W9Bl . The straight line is p,, = S/-^. 

to write GC densities ab initio in terms of the (fixed) peri- 
centers rp of unique orbits in the adopted potentials. Cluster 
relaxation times and mass-loss rates /igv then follow as func- 
tions of rp as well. Finally, specific initial mass, space, and 
velocity (or orbital eccentricity) distributions are chosen for 
entire GC systems, so that at all later times it is known what 
the dynamically evolved dN/ d log M is for globulars with any 
single rp-, how many clusters with a given rp survive; and what 
the distributions of rp and all dependent cluster properties are 
at any instantaneous position r^^. 

In this approach, if the GCMF began with a power-law rise 
towards low masses and its current peak is due entirely to clus- 
ter disruption, then a dependence of Mjo on rp is expected in 
general, because the densities of tidally limited GCs decrease 
with increasing rp. Thus, models along these lines that as- 
sume the orbit distribution of a GC system to be the same 
at all radii in a galaxy (i.e., that the time average of the ra- 
tio r^c/rp is independent of position) have typically had diffi- 
culty in accounting for the observed weak or non-correlation 
between Mto and present r^^ in large galaxies. This is partic- 
ularly a problem if it is assumed that the initial GCMF was a 
pure power l aw, with the same index at arbitraril y high masses 
as low (e.g., lBaumgardllll998t IVesperinill2001h . I t is poten- 
tially less of a concern if dN /d logM started as a ISchechteJ 
function with an exponential cut-off at masses M > Mc, as 
we have assumed, since then the existence of a strict upper 
bound Mjo < Mc 02.21 1 means that the dependence of an 
evaporation-evolved Mto on rp and rgc must saturate for small 
enough galactocentric radii (high enough GC densities). Even 
so, the "scale-free" models of FZOl, in which ~ 10^ Mq 
and all GCs in a Milky Way-like galaxy potential have the 
same time-averaged r^^/rp, still predict a gradient in Mto ver- 
sus r^c. that is stronger than observed. 

iFZOll showed that, if they left all of their other assumptions 
unchanged, then a dependence of GCMF peak mass on rgc 
could be effectively erased by an appropriately varying radial 
velocity anisotropy in the initial GC system. Thus, in their 
"Eddington" models the eccentricity of a typical cluster or- 



bit increases with galactocentric distance (the time average of 
fgc/rp increases with radius), such that globulars spread over 
a larger range of current rgc can have more similar rp and asso- 
ciated Mto- However, the initial velocity-anisotropy gradient 
required to fit the Milky Way GCMF data specifically is only 
marginally c onsistent with the observed kinematic s of the GC 
system ( e.g..lDinescu. Girard. & van Altenal[T999h .'^ Subse- 
quently, IVesperini et al.l (|2003|) constructed broadly similar 
models for the GCMF of the Virgo elliptical M87 and con- 
cluded that there, too, a variable radial velocity anisotropy is 
required to match the observed Mto versus rgc; but the model 
anisotropy profile in this case is clearly inconsistent with the 
true velocity distribution of the G C system, which is observed 
to be isotropic out to large rgc (iRomanowskv & Kochanekl 
'2Ml lC6tietani200l . 

These results certainly suggest that some element is lacking 
in rgc-oriented GCMF models developed as outlined above. 
But they do not mean that the fault lies with the main hypoth- 
esis, that the difference between the mass functions of young 
clusters and old GCs is due to the effects of slow, relaxation- 
driven disruption in the latter case. Any conclusions about 
velocity anisotropy depend on the totality of steps taken to 
connect the densities and positions of clusters; and it is possi- 
ble that reasonable changes to one or more of these ancillary 
assumptions could make the models compatible with the ob- 
served kinematics of GCs in both the Milky Way and M87, 
without abandoning a basic physical picture of evaporation- 
dominated GCMF evolution that is otherwise quite success- 
ful. 

One issue is that previous models have always specified 
evaporation rates a priori as functions of cluster density (or or- 
bital pericenter), usually normalizing /^ev so that fdis/frh — 20- 
40 as in standard treatments of t wo-b ody rela xation. How- 
ever, following our discussion in ^3.1l and ^3.21 it would seem 
worthwhile to investigate these models with /igv increased at 
fixed pk or rp to allow fdis/frh ~ 10 (if /3 ~ 2 for the low-mass 
p ower-la w part of the initial GCMF). 

iFZOlland lVesperini et al.l (12003) both consider velocity dis- 
tributions parametrized by a galactocentric anisotropy radius, 
Ra, inside of which a cluster system is essentially isotropic 
and beyond which it is increasingly dominated by radial or- 
bits. In these terms, the difficulty with the published models 
is that, to reproduce the observed insensitivity of Mto to rgc 
given standard normalizations of /iev, they require values of 
Ra that are smaller than allowed by observations (especially 
for M87). Increasing Ra to more realistic values while keep- 
ing the normalization of /igv fixed leads to a stronger gra- 
dient in Mto: the orbits of GCs at small rgc < Ra remain 
closely isotropic and the typical rp and Mto are essentially 
unchanged, while at large galactocentric distances the clus- 
ter orbits are on average less radial than before, with larger 
rp, lower densities, and lower evolved Mto for a given rgc. 
This effect is illustrated, for example, in Figure 9 of FZOll 
However, it can be compensated at least in part by increas- 
ing /iev by a common factor for all GCs, with the new, larger 
R a fixed, if t he initial mass function is assumed to have been 
a \Schechten function rather than a pure power law extend- 
ing to arbitrarily high masses. A faster evaporation rate 
will then lead to a (roughly) proportionate increase in the 

" The fact that clusters on radial orbits are preferentially disrupted lessens 
any inconsistency between the radial anisotropy required in the initial veloc- 
ity distribution and observational constraints on the present velocity distribu- 
tion. 
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evolved GCMF peak mass for GCs with relatively low den- 
sities, i.e., those at large r^c and r^; but the increase in Mto 
will be smaller, and eventually even negligible, for higher- 
density clusters at progressively smaller rgc — again because 

1/2 

Mto grows less than linearly with fi^v oc p,/ when there is 
an upper limit Mjo < Mc due to an exponential cut-off in the 
initial dN/d log Mq. Thus, the qualitative effect of increasing 
the normalization of /iev in models with radially varying GC 
velocity anisotropy is to weaken the amount of radial-orbit 
bias required to fit an observed M m vers us r^^. 

Another point, emphasized by FZOl', has to do with the 
standard starting assumption that GCs orbit in galaxies that 
are perfectly static and spherical. In reality, galaxies grow 
hierarchically. In this case, even if the values of /Xev are not 
changed, much of the burden for the weakening or erasing of 
any initial gradients in Mjo versus may be transferred from 
velocity anisotropy to the time-dependent evolution of the 
galaxies themselves. Violent relaxation, major mergers, and 
smaller accretion events all work to move clusters between 
different parts of galaxies and between different progenitors, 
scrambling and combining any number of pericenter-density- 
Mto relations. Any position dependences in the GC ph dis- 
tribution and in Mto itself for the final galaxy are therefore 
bound to be weaker, more scattered, and more difficult to re- 
late accurately to a cluster velocity distribution than in the 
case of a monolithic, non-evolving potential. Allowing for 
a non-spherical galaxy potential would have qualitatively the 
same effect, because in this case every cluster explores a range 
of pericenters and different maximum tidal fields on each of 
its orbits. 

In this situation, it may be important to ask how evapora- 
tion rates can still be approximately constant in time — so that 
cluster masses still decrease approximately linearly with t as 
our models assume — if the tidal field around any given GC 
changes significantly over time. Thus, consider first a sys- 
tem of GCs in a single, static galaxy potential. The mass- 
evolution curve for each cluster is approximately a straight 
line, M(f ) ~ Mq - pieyt, with /igv depending on some measure 
of internal density, which may be pj/^, p]^'^, or Y^J'^ . The av- 
erage mass-evolution curve for the entire system of clusters 
is also approximately linear, {M(t)) ~ (Mq) - {pev)t. If now 
a merger or other event rearranges the clusters in the galaxy, 
then after the event the mass-loss rates of some clusters will 
be higher than before and the rates of other clusters will be 
lower than before. However, if the mean density of the galaxy 
as a whole is roughly the same after the event as before, then 
so too will be the average of the GC densities, because of 

1 /2 

tidal limitation. The average (/iev) oc (p,/ ) (say) will differ 
even less between the pre- and post-merger systems. Thus, al- 
though using instantaneous densities to estimate the past /igv 
of individual clusters may err on the high side for some clus- 
ters and on the low side for others, these errors will average 
away to a small or even zero net bias. The approximation 
/Xev — constant in time in our GCMF models will then still be 
valid in the mean, and the average (M(f)) dependence of suf- 
ficiently large numbers of clusters will remain roughly linear. 

This type of scenario might be expected to pertain at least 
to galaxies that evolve on the fundamental plane, since this 
entails a connection between the total (baryonic plus dark) 
masses and circular speeds of galaxies, of the form Mgai oc 
or Mgai oc V^. By the virial theorem, the average densities 
scale as pgai oc V^/M^, and thus pgai oc M^^, or pgai oc M~y^. 



Insofar as (pi,) oc (p,) oc pgai for the GCs, the system-wide av- 

1 /2 

erage (/iev) cx (p,/ ) should therefore not change drastically 
even after a major merger between two fundamental-plane 
galaxies; at most, the ratio of final to initial (/iev) will be 
roughly of order the -1/4 power of the ratio of final to ini- 
tial Mgai. Note that thi s line of reasoning is closely related to 
that applied by iFZOll to explain the small observed galaxy- 
to-galaxy differences in the average turnover masses of entire 
GC systems (although non-zero differences do e xist, and can 
be accomod ated in these sorts of arguments; see iJordan et al.l 
l2006ll2007h . 

A full exploration of questions such as these, about the wide 
range of ingredients in current GC-plus-galaxy models, will 
most likely require large A^-body simulations set in a realistic, 
cold dark matter cosmology. Until these can be carried out, 
it is our view that the kinematics of globular cluster systems 
cannot be used as decisive side constraints on theories for the 
GCMF 

4. CONCLUSIONS 

We have shown that the mass function dN/d log M of glob- 
ular clusters in the Milky Way depends significantly on clus- 
ter half-mass density, ph, with the peak or turnover mass Mto 
increasing and the width of the distribution decreasing as ph 
increases. This behavior is expected if the GCMF initially 
rose towards masses below the present turnover scale — as the 
mass functions of young cluster systems like that in the An- 
tennae galaxies do — and has evolved to its current shape via 
the slow depletion of low-mass clusters over Gyr timescales, 
primarily through relaxation-driven evaporation. The fact that 
Mto increases with cluster density favors evaporation over 
external gravitational shocks as the primary mechanism of 
low-mass cluster disruption, since the mass-loss rates asso- 
ciated with shocks depend inversely on cluster density and 
directly on cluster mass. Our results therefore add to previ- 
ous arguments supporting an interpretation of the GCMF in 
terms of evaporation-dominated evolution, based on the fact 
that d N/d log M scales a s M'"'^ with /3 ~ in the low-mass 
Hmit (.Fall & Zhangll2001h . 

The observed GCMF as a function of pi, is fitted well 
by simple models in which the initial distribution was 
a ISchechteri function, dN/d log Mq oc m]^^ exp (-Mq/M,.) 
with (5 = 2 and M^. ~ 10* Mq assumed, and in which clusters 
have been losing mass for a Hubble time at roughly steady 
rates that can be estimated from their current half-mass den- 

1 /2 

sities as /iev k Ph ■ We have shown that, although this pre- 
scription is approximate, it captures the main physical depen- 
dence of relaxation-driven evaporation. In particular, it leads 
to model GCMFs that are entirely consistent with those re- 
sulting from alternative characterizations of evaporation rates 
in terms of cluster tidal densities p, or mean surface densities 
( ^3.2l i. The normalization of /iev at a given ph (or p,, or S,) 
required to fit the GCMF implies total cluster lifetimes that 
are within range of the lifetimes typically obtained in theoret- 
ical studies of two-body relaxation, although our values may 
be slightly shorter than the theoretical ones if the low-mass, 
power-law part of the initial cluster mass function was as steep 
as we have assumed. 

Taking clusters in various bins of central concentration c 
and Galactocentric radius r„c and using their (individual) ob- 
served densities as direct input to our models yields dynam- 
ically evolved GCMFs as functions of c and r^^ that agree 
well with all data. This again indicates that the most fun- 
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damental physical dependence in the GCMF is that on clus- 
ter density. Moreover, our models for dN /d log M versus 
obtained in this way are consistent in particular with the 
well-known insensitivity of the GCMF peak mass to Galac- 
tocentric position. This is seen to follow from a significant 
variation of Mjo with pf, (or p,, or E,) — due in our analysis 
to evaporation-dominated cluster disruption — combined with 
substantial scatter in the GC densities at any Galactocentric 
position. 

We have not invoked an anisotropic GC velocity distribu- 
tion to explain the observed weak variation of Mjo with r^^', 
indeed, we have made no predictions or assumptions what- 
soever about velocity anisotropy. We have emphasized that, 
when velocity anisotropy enters other long-term dynamical- 
evolution models for the GCMF, it is only in conjunction with 
several additional, interrelated assumptions made as part of 
larger efforts to derive theoretical density-rgc relations for 
GCs — which we have not attempted to do here. The appar- 
ent need in some current models for a strong bias towards 
high-eccentricity cluster orbits to explain the near-constancy 
of Mjo versus rgc might well be avoided by changing one or 
more ancillary assumptions in the models, without having to 
discard the underlying idea that the peak and low-mass shape 
of the GCMF are the result of relaxation-driven cluster dis- 



ruption. 

It clearly will be of interest to test and refine the main 
ideas in this paper through modeling of the GCMFs in other 
galaxies. For the time being at least, doing so will re- 
quire the estimation of approximate mass-loss rates using 
cluster half-mass densities rather than tidal quantities, sim- 
ply because GC half-light radii can be measured accurately 
in many systems beyond the Local Group, whereas tidal 
radii ar e much more model-dependent and difficult to ob- 
serve. IChandar. Fall. & M cLaughlinI (l2007h have recently 
shown that the peak mass of the GCMF in the Sombrero 
galaxy (Ml 04) increases with p/, in a way that is reasonably 
well described by sums of evolved lSchechted (1 1 976b functions 
as in the models presented in this paper. It should be rela- 
tively straightforward to pursue similar studies in other nearby 
galaxies. 
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